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Surface excitations in electron spectroscopy. 
Part I: dielectric formalism and Monte Carlo 
algorithm 

F. Salvat-Pujol'' and W. S. M. Werner 

The theory describing energy losses of charged non-relativistic projectiles crossing a planar interface is derived on the basis of the 
Maxwell equations, outlining the physical assumptions of the model in great detail. The employed approach is very general in that 
various common models for surface excitations (such as the specular reflection model) can be obtained by an appropriate choice 
of parameter values. The dynamics of charged projectiles near surfaces is examined by calculations of the induced surface charge 
and the depth- and direction-dependent differential inelastic inverse mean free path (DIIMFP) and stopping power. The effect of 
several simplifications frequently encountered in the literature is investigated: differences of up to 100% are found in heights, 
widths, and positions of peaks in the DIIMFP. The presented model is implemented in a Monte Carlo algorithm for the simulation 
of the electron transport relevant for surface electron spectroscopy. Simulated reflection electron energy loss spectra are in good 
agreement with experiment on an absolute scale. Copyright © 2012 John Wiley & Sons, Ltd. 
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Introduction 

A quantitative understanding of any spectroscopy technique 
based on the analysis of reflected, transmitted, or emitted charged 
projectiles from solid surfaces relies on an accurate description of 
the inelastic interaction of the projectiles with the target. In the 
bulk of the solid, the stopping of the projectile can be accurately 
described within the semiclassical dielectric formalism for infinite 
media: the presence of the charged projectile perturbs the equilib- 
rium charge distribution of the solid, which becomes polarized 
and induces an electromagnetic (e.m.) field that acts as a stopping 
force on the projectile. Within the semiclassical approximation, the 
Fourier spectrum of the field set up by the projectile is interpreted 
as a distribution of energy losses and momentum transfers at 
individual, discrete inelastic interactions. Thus, the semiclassical 
approximation serves as the link between the classical dielectric 
description and the quantal nature of inelastic interactions. The 
presence of an interface, e.g. solid-vacuum or between two solids, 
imposes strict boundary conditions on the total (and hence also 
on the induced) e.m. field. Thus, the stopping of the projectile close 
to the interface is substantially different than in the bulk of the solid: 
an additional dependency is found on the distance to the interface 
(at either side) and on the surface crossing direction. 

A collective response of the electrons in the surface region of 
the solid due to the passage of charged projectiles through the 
interface was predicted in the 1950s by Ritchie'^^ and confirmed 
experimentally by Powell and Swan.'^'^' A number of models have 
been developed in the last decades to understand the collective 
(and single-particle) response.'^^^' Different approaches and 
approximations are adopted: some models use a classical electro- 
dynamics framework whereas others use many-body quantum 
theory, some assume a simplified dielectric response of the solid, 
some are valid for particular trajectories of the projectile (e.g. 



perpendicular or parallel to the surface). Simplifying mathematical 
assumptions are often made in order to highlight the relevant 
physics, to obtain more treatable expressions, and to keep the 
computation time within reasonable limits. Unfortunately, the effect 
of the different approximations on the resulting physical description 
of the stopping of the projectile has not been scrutinized. 

It is the purpose of the present work to present a self-contained 
description of the dynamics of charged projectiles near a planar 
interface within the semiclassical dielectric formalism, taking 
special care to avoid further approximations such as those 
addressed above. A common framework is provided which allows 
one to compare a subset of models in the literature and to expose 
the effects of their implied approximations. The calculation 
procedure presented here has been implemented in a Monte Carlo 
algorithm for the simulation of reflection-electron-energy-loss 
spectra. Finally, the quality of the model will be assessed by means 
of comparisons of simulated reflection electron energy loss spectra 
(REELS) with experimental REELS on an absolute scale. 

This work is structured as follows: In section "Theory", the dielec- 
tric formalism for the stopping of charged projectiles in an infinite 
medium is briefly reviewed and extended to take into account 
the presence of a planar interface to another medium. An expres- 
sion for the differential inelastic inverse mean free path (DIIMFP) 
of a charged projectile in the vicinity of a planar interface is derived. 
In section "Results", the dependence of the derived DIIMFP on the 
velocity of the projectile, on the distance at each side of the 
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i) Semi-infinite geometry 
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Figure 1 . (i) Considered geometry: a non-relativistic projectile of charge Zoe and velocity v is at a distance d at either side of the interface between two 
media tJs\,o)) and fb(q,oj) moving with a polar angle ;/. (ii) and (iii): Schematic representation of the method of image charges for dielectrics. By an 
appropriate choice of the contributions p,, p2, and p3 (unity or zero) different models found in the literature can be reproduced, such as the surface- 
reflection model, for which Pi =P2 = P3 = 1 [see text and Eqn (8)]. 



Interface, and on the surface crossing direction is examined. 
Calculated DIIMFPs are compared with those given by another 
model from the literature, exposing the effects of its implied 
approximations. Next, a IVlonte Carlo algorithm for the simulation 
of electron reflection experiments is described. Simulated spectra 
are compared in absolute units with experimental spectra. Finally, 
in section "Conclusions", the presented results are summarized, 
and conclusions are drawn. 



Theory 

We consider the geometry depicted in Fig. li. A projectile of charge 
ZqC moves near the interface between two semi-infinite media with 
dielectric functions ejs\fii) and fb(t|/'j)/ where the subscripts a and 
b imply above and below the interface, respectively.^ For the formu- 
lae to apply for electrons, Zq must be set to -1 . We use a system of 
Cartesian coordinates where the surface is the plane z = 0 and the 
unit vector along the z axis, z, points towards medium a. The 
charged projectile moves in a uniform motion with velocity v 
making an angle r; with the positive z axis. We denote the z coordi- 
nate of the projectile's position as d. There is mirror symmetry with 
respect to the plane spanned by v and the z axis, and, furthermore, 
in the case of normal escape and normal incidence (v parallel and 
opposed to z, respectively), all physical quantities we derive must 
be invariant under an arbitrary rotation around the z axis. We will 
use the Gaussian system of units [l/(47ifo) = 1, where is the 
vacuum permittivity, and ;(o/(47r) = l, where /io is the vacuum 
permeability] throughout the calculation, unless stated otherwise. 

First, the method of image charges will be used to solve 
Maxwell's equations in order to obtain the total e.m. field set 
up by the projectile with the appropriate boundary conditions. 
Next, the e.m. field induced at the position of the projectile will 
be determined. Finally, an expression for the average energy loss 
per unit path length (the so-called electronic stopping power) 
will be derived as a function of the distance to the interface, d, 
the surface crossing angle, r/, and the speed of the projectile, v. 
The connection between this classical description and the quan- 
tized nature of the inelastic interaction of the charged projectile 
in the solid will be made at the final stage of the calculation by 

^Although in practice the medium above the surface will often be taken as 
vacuum, with «o(q,tu) = 1 , we keep an arbitrary medium a to retain generality. 



means of the semiclassical approximation: and ftco will be con- 
sidered as the momentum and energy transfers from the projec- 
tile to the solid in an individual inelastic interaction. Here, q and 
w are the Fourier conjugate variables of the position and time 
variables r and f, respectively. By strict analogy with the quantal 
calculation, one may derive not only the average energy loss 
per unit path length, but also the DIIMFP, a quantity proportional 
to the probability of undergoing a discrete energy loss between 
1/1/ and l/l/-tdl/l/ per unit path length. 

We will work in the Coulomb gauge, where the vector poten- 
tial A(r,t) is transverse: V A(r,f) = 0. We recall that the electric 
field in terms of the scalar potential ip(t,t) and the vector potential 
A(r,f) is given by^^"^' 



-V</p(r,r) 



1 aA(r, t) 

c dt 



(1) 



where c is the speed of light in vacuum. As we are concerned 
with energy losses of relatively slow projectiles (with speeds 
much less than c), the transverse components of the e.m. field 
can be ignored. Consequently, we take A = 0, so that (B is thus 
automatically set to zero): 



f(r,t) = -V^(r,t). 



(2) 



The effect of this approximation is twofold: (i) Maxwell's 
equations are formally reduced to the electrostatic case even 
though the position of the charged projectile and, hence, the 
induced field vary with time and (ii) the theory does not describe 
the effect of 'transverse interactions', which contribute to the 
energy loss of relativistic projectiles and cause effects such as 
transition radiation and Cerenkov radiation. In the case of 
electrons, this implies that the present calculation is valid for 
kinetic energies up to a few tens of keV. Retardation effects due 
to the finite propagation speed of the e.m. field are not accounted 
for explicitly. 

We assume that the perturbation caused by the passage of 
the projectile through the media is small enough so that the 
media respond in a linear fashion. For homogeneous media, 
this is equivalent to assuming that the electric displacement 
field X>(r, f) and the electric field are related via the linear 
relationship 



I>(r, f) 



df £(r - r , t - f )f (r ,f ' 



(3) 
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In Fourier space, this expression becomes 
D(q,a;) = e(q,co)f(q,w), 



(4) 



where e(q,co) is the dielectric function of the medium. The linear- 
response approximation is formally equivalent to a first-order Born 
approximation in perturbation theory, which is valid for projectile 
energies of about 20 times the characteristic energies of the 
active target electrons (ionization energies in the case of inner-shell 
excitations).'^''' Other estimations in the case of inelastic scattering 
in solids'^^' suggest a lower bound of valid kinetic energies in the 
order of 7Ef, where £p is the Fermi energy of the solid. Thus, in 
the case of electrons, the present calculation is strictly valid for 
kinetic energies between ~ 1 00 eV and ~ 1 0 keV. 

We recall that the dielectric function gives the response of 
the solid to the {q,oj)-Fourier component of an external e.m. field. If 
the external perturbation is a plane wave with frequency to 
propagating along a given direction, e.g. the electric field carried 
by incoming photons in the case of x-ray or UV spectroscopy, the 
modulus of the electric field is independent of the spatial coordinates 
and, hence, the response of the solid is characterized by a dielectric 
function that depends only on the frequency of the wave, f (to). If 
the perturbation is due to a moving charged projectile, the electrons 
of the solid around the projectile interact with it via the Coulomb 
potential, which evidently depends on the spatial coordinates. Thus, 
the polarization of the material is non-local and is therefore in 
general a non-trivial function of both r and r, or q and co in Fourier 
space. Therefore, strictly speaking, in order to describe the dielectric 
response of a solid to the passage of a charged projectile, a dielectric 
function of the form e(q,oj) is needed. This kind of dielectric function 
is known in the literature as spatially dispersive, owing to the 
dependency on q, namely on r in real space. 

Finally, we recall that the boundary conditions of the electric field 
and of the electric displacement field at the interface read'^'^^' 



limz^o+f||(»',t) 
(r,r) 



lirriz^o- f||(r,f) 



(5) 



where the subscripts || and _L indicate, respectively, the compo- 
nents parallel and perpendicular to the surface. These boundary 
conditions must hold at all times at any point on the surface. They 
also apply for time-varying fields.'^^' 

Method of image charges 

Maxwell's equations for electrostatics with a given set of 
boundary conditions can be formally solved by means of the 
Green's function formalism.'^*'^^' The potential is given as an 
integral relationship within a closed volume involving a kernel of 
the form G{t,t) = 1/|r — r'| + F(r, r'), where r is the observation 
point and t is the source point. The first term of G[t,t) is the 
potential at r due to a point charge at t . The second term, f (r, r ), 
is seen to satisfy the Laplace equation'^^' {V ^</?(r,t) = 0, potential 
in a region of space with no charges) and can therefore be 
understood as the potential of a system of charges external to the 
integration volume, chosen in such a way that, together with a 
source at r', the required boundary conditions are satisfied. The 
method of image charges for dielectrics,'^^' also known in the 
literature as the method of extended pseudomedia,'^^'^'' is a prac- 
tical approach which is formally equivalent to determining F(r, r ), 
with the benefit of a conceptually simpler development. It consists 
in finding an alternative representation of the problem in terms of 
fictitious charge distributions chosen in such a way that the electric 



field (i) can be trivially found and (ii) verifies the same differential 
equation (Maxwell's equations) and the same boundary conditions 
as in the original problem. Such an alternative representation is 
therefore equivalent to the original problem for all purposes. 

We consider two infinite pseudomedia a and b with bulk 
dielectric functions fa(q,(o) and fi,(q,co), respectively (see Fig. lii 
and iii). In order to emulate the effect of the interface in the real 
problem, fictitious charge distributions will be added to each 
pseudomedium. Next, the electric field in each pseudomedium, 
Ea{t,t) and Eb{tA), will be determined. The electric field in the 
semi-infinite geometry will then be given by 



f) for d > 0, 
r) for d < 0. 



(6) 



At this point, the electric field will still contain the fictitious charge 
distributions explicitly. Next, the boundary conditions (5) will be ap- 
plied and the fictitious charges shall be completely determined in 
terms of other parameters and disappear from the equations, thus 
serving as a convenient means of solving Maxwell's equations sat- 
isfying the corresponding boundary conditions at the planar inter- 
face. Each infinite pseudomedium exhibits translational invariance: 
its dielectric response can be described in terms of the bulk dielec- 
tric function. Thus, the method of image charges allows us to de- 
scribe surface effects using bulk properties of the media (I). We 
point out for completeness that an alternative approach consists 
in determining the surface dielectric response explicitly.'^°'^^' 

The fictitious charge distributions at each pseudomedium will 
now be specified. The question arises as to whether the electric 
field is uniquely defined, that is, whether different choices of 
fictitious charges yield the same electric field. We recall that the 
electrostatic problem is uniquely determined only for Dirichlet or 
Neumann boundary conditions'^^' (specifying the potential or the 
normal derivative of the potential, respectively, everywhere on a 
closed surface). The boundary conditions, Eqn (5), are neither 
Dirichlet nor Neumann boundary conditions: only continuity of 
the parallel component of the electric field and of the normal 
component of the displacement through the surface is required, 
not a specific value. We therefore assume that the electric 
field which satisfies the required conditions is not unique. Other 
authors make similar statements: Chan et a/.'"' state that 'different 
choices of fictitious charge distributions correspond to different 
boundary conditions' and, similarly, Flores et o/.'^^' mention that 
'there exists a precise relationship between boundary conditions, 
physical models of surface scattering, and corresponding surface 
pseudostimuli.' Thus, it must be observed that different choices of 
fictitious charges yield different descriptions of the inelastic 
scattering of charged projectiles near a surface. See the note after 
Eqn (8) for a further clarification of this point. 

Recall that the velocity of the projectile is given by y = (v„Vy,v^. 
The planar interface z=0 suggests a candidate fictitious charge at 
each pseudomedium: the specular image charge of the projectile, 
which moves with a velocity v = {vx, Vy, —Vz). Furthermore, as we 
will see below, the use of spatially dispersive dielectric functions 
implies that a fictitious surface charge ps(t\\,t) must be added at 
the plane z=0.'^^'^^' Thus, we express the charge distributions in 
pseudomedia a and b as 



p„(r, f) = Zoe&(t - Mt) + p,Zoe&{t - vt)+ p,(r, t) 
Pb(r, t) = P2Zoe8{t - vt) + P3Zoe<5(r - vt) - p^{r, t) 



(7) 



respectively. Notice the opposite sign of the surface charge in 
pseudomedium b in Eqn (7). This sign is necessary; otherwise, the 
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normal component of the electric displacement field cannot be 
continuous through the surface.'^^' The scalar parameters pi, P2, 
and p3 have been introduced to conveniently reproduce different 
models found in the literature: 



(Pl,P2,P3) 



(1,1,1) Ref.'*'2^'^''(specular reflection model, SRM) 

(0,1,0) Ref.i«'2i.23 

(8) 

If the calculation is carried out including the transverse part of the 
electric field, one finds that the boundary conditions for the trans- 
verse part of X>(r, t) can only be satisfied if the image charge is se- 
lected, e.g. (pi, p2. Pa) = (1,1,1). It is then seen that this point image 
charge accounts for a bulk polarization current moving away from 
the interface as the projectile moves towards it in the medium side. 
It is implicitly understood that the electrons of the solid bounce 
back towards the inside of the medium at the interface as they 
are repelled by the projectile (assuming an electron). This picture 
corresponds to the specular reflection model.™ In the present cal- 
culation, where we neglect the transverse part of f (r, t), the bound- 
ary conditions are not restrictive enough: they allow for the variety 
of image charges (and therefore of surface-scattering models) dis- 
cussed above. 

Figures lii and liii show p^ir, t) and pj,(r, t), respectively, for a 
case where the charge is in medium a. The Fourier transforms 
[Eqns (A1) and (A3)] of the charge distributions read 



Po(q. cj) = ^ [S{co - q-v) + pi(5(a) - q-v)] + p,(q, w) 

2.K 

Pb{<\, co) = ^ [P2<5(co - q-v) + p3<5(co - q-v)] - p,(q, w). 
zn 



(9) 



Taking the Fourier transform of the boundary conditions (5) and us- 
ing Eqns (13) and (14), we have that the boundary conditions in 
terms of po(q, w) and pb(q, co) are given by 



J — c 



c?^ + (7|e<,(q,w) 



(15) 



and 



dc?z 



f 

J — c 



respectively. To derive the last two expressions, we have used 
the fact that the boundary conditions must hold everywhere 
on the surface at all times (Vx,y, f). Notice that the dielectric 
functions at each side of the interface depend on q and co. For 
the equalities to hold for all qy and m, the charge densities must 
have a continuous dependency on q and co. A collection of point 
charge distributions of the form p{q,co) ~ (5(co — q - v) is therefore 
not enough to satisfy these two equalities for dispersive dielec- 
tric functions. An extended charge distribution with dependency 
on q|| and co (a surface charge distribution) is required, as we 
anticipated above. 

The boundary conditions are now imposed. We start by 
inserting po(q, co) and pf,(q, co), Eqn (9), into Eqn (1 5). The fictitious 
surface charge distribution which renders the two infinite 
pseudomedia consistent with the semi-infinite geometry can be 
immediately isolated: 









} 






1 

Ca(q||,Kz,tO 







We now express the boundary conditions in terms of pa(q, co) where in the numerator, q should be understood as (qy, /CJ. The 
and pf,(q, co). First, the electric field and the displacement field will integral over can be readily carried out using 
be expressed in terms of Pa(q, cu) and Pb(q, co). Poisson's equation 
reads'^='^« 



V-I>(r, t) = 47ip(r, f) 



(10) 



.(co-q-v)=^. /Cz-'^^ 



(18) 



where p(r, f) is the distribution of external charges. In Fourier space, 
this equation becomes 



iq-X>(q, co) = 47ip(q, co). 



(11) 



where p{r, f) is the distribution of free charge. Using Eqn (4) and the 
Fourier transform of Eqn (2), 



f(q,&j) = -iq'j=(q,&j 

we can recast Eqn (1 1) as 

q 1 



£(q,co) 



-1471 



c?2e(q,co) 



p(q,«j 



and write 



D(q,co) 



q / N 

-1471 — p(q.co) 



(12) 



(13) 



(14) 



Defining 



and 



1 / CO — qii-vi 
,5(co-q-v) =— (51 /<z+ " 



Vz 



Pi 



1 



P3 



(19) 



(20) 



Pi 



ct,{f\pK,ioj Ca{<\pk2,oi^ ct,[q\\,-k^.oi^ c„ (q|| , -/c^, to) 
1 



1 

q» + «l 



1 



1 



£0 qih'i'z.co fb qih^z.co 



(21) 
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we can recast Ps(q,co) as 



Ps(q-«) 



271 Cj2^2 + (cO-q||-V||)2 



f(qii,co) 



(22) 



This expression coincides witli Eqn (3) in RefJ^*' witliin notation 
changes. The second boundary condition, Eqn (1 6), yields a closure 
relationship between p,, p2, and which restricts the value of the 
fictitious point image charges: 



Pi + P2 - P3 = 1 • 



(23) 



All literature models listed in Eqn (8) satisfy this relationship. 

Once ps(q, oj) has been obtained, the total electric field, Eqn (1 3), 
is completely determined as 



£(q,(u) = -l47r^ 



" Q(d) /Zpe 
c„(q,a)) I 271 

, w) I In 



3{aj - q ■ v) + pi r5((o - q ■ v)] + ^^(q, co) 



+ "I ^ IP2&{oj - q ■ v) + P3<5((;j - q ■ v)] - p,(q, co) 

Cfclq,'") 271 



(24) 



which can be recast more conveniently as 



£{q,co) 



0(d) , P20(-d) 



+ —d(oj-q-v 
In 



eo(q,co) eb(q,(u) 
Pi0(d) , P30(-d)" 



ea(q,wj eb(q,coj 
0(d) 0(-d)- 



eo(q,a)) e6(q,co) 



(25) 



where ©(d) is the Heaviside step function, which is 1 if d> Oand 0 if 
d<0. The induced electric field, find(q, «), is obtained by subtract- 
ing the field that the projectile would create in vacuum, fvac(q, co), 
from the total electric field. From Eqn (13), we have that 



£vac(q,co) = -i47i^g,5(a;-q.v) 



(26) 



Thus, the induced electric field is given by 



£ind(q,a') = -1471 



I 271 



8{o) — q-v) 



0(d) , p20(-d) 



271 



6{w — q-v) 



(a(q,(o) eb(q,w) 
p,0(d) , P30(-d)l 



fa(q,<») eb(q,«) 



Q(rf) Q(-rf) 

eo(q,w) efa(q,oj) 



(27) 



Finally, using Eqn (12), we have that the induced potential is 
given by 



271 



0(d) 



P2@(-d) 



Zoe 

' 271 



S(co — q-v) 



eo(q,cu) eb(q.(u) 
p,0(d) , p30(-d)' 



-P5(q,cj) 



0(d) 



eo(q,fj) f(,(q,oj) 



0(-d) 



fo(q,co) f6(q,co) 



(28) 



The induced electric field can be interpreted as the superposition of 
a field created by induced charges in the bulk of the media [first 
two terms in Eqn (27)] and a field created by an induced surface 
charge [third term in Eqn (27)]. Thus, the fictitious surface charge 
Ps(r, t), which has so far served as a trick to satisfy the boundary 
conditions, is found to coincide with the induced surface charge. 

Induced surface charge 

Before carrying on with the calculation of the stopping power, it 
is instructive to examine in some detail the induced surface charge 
Ps(r, t) imposed by the boundary conditions, since all surface-related 
quantities will be obtained as integral transformations of p^ir, f) 
[see, e.g. Eqns (25) to (28)]. Thus, knowledge of the dependency 
of ps(r, t) on d, r/, and v gives valuable insight into the conditions 
for which surface effects are to be expected. 
The inverse Fourier transform of Eqn (22) reads 



< rdcof dqii e'(''« ''ir™0 



2^^CO-q||-V||)2' 

(29) 

where we have carried out the trivial integration over q^ and where 
the term sign(d) has been explicitly introduced here to automatically 
select the appropriate sign of ps(q, co) depending on whether the 
projectile is above or below the surface, see Fig. 1 and Eqn (7). With 
the change of variables £3 = co — qii - vm, we have 



1 

^5(''>f) = ^^oekz|5(z)sign(d) 



dqii 



dco e 



-i«t''(qii'«+qirvii 



qlvj + 



(30) 



Causality arguments imply that l/c(q,oj) is analytic in the upper half- 
plane of complex co.'^^' Thus, f(q||,co) has all poles in the lower half- 
plane. For f < 0, we can therefore perform the integration over co by 
closing the integration contour as indicated in Fig. 2 and using Eqns 
(B3) and (B8) to obtain 

p,(r,f) = -^Zoe,5(z)0(-t)sign(d) /"dc?,, e-""!"^! 
27r Jo 

X y^"^^'^d,^'Re[e'''ir('ii-''iiOf(q||,q||.V|| + ic?|||u,|)], 

(31) 

where q^ =c?||(cos(^',sin0',O). 

In the case f > 0, we consider Eqn (30) and use the second line 
of lemma (Bl) and the property e(q,co) = e"(— q, — co) to obtain 
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Figure 2. Integration contour used to carry out the integral over co in 
Eqn (30). 
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(32) 
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Figure 3 displays the surface polarization charge p^{t, t), 
Eqns (31) and (32), induced by an outgoing 500 eV electron that 
moves from Al (d<0) to vacuum (c/>0) perpendicular to the 
interface (77 = 0) for different values of the distance to the inter- 
face, d. The dielectric-function model used to produce this figure 
is given below, Eqn (66). The figure should be read by columns in 
order to follow the time evolution of ps(r, t) at different stages of 
the trajectory of the electron. The black arrows with the label v 
represent the direction of motion of the electron. When the 
electron approaches the interface from the Al side (d<0), the 
induced surface charge ps(r, r) is negative, and therefore a 
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repulsive force is induced on the electron. This force, 
schematically represented in the figure by the arrows labeled 
Find.s/ opposes the direction of motion of the electron: the 
electron is decelerated. Similarly, when the electron is at the 
vacuum side (d> 0), psir, f) is positive, and therefore an attractive 
force is induced on the electron. This attractive force opposes 
the direction of motion of the electron: the electron is also 
decelerated at the vacuum side. Thus, the contribution of the 
induced surface charge to the total induced force on an outgoing 
projectile moving perpendicularly to the surface results in a 
deceleration at both sides of the interface. 

In Fig. 4, we consider the opposite trajectory: an incoming 
electron with an energy of 500 eV crosses the interface from 
the vacuum side (rf>0) to the Al side (d<0) in a direction 
perpendicular to the surface (r/=180°). As the electron 
approaches the surface from the vacuum side (d>0), a positive 
charge is induced at the surface, which creates an attractive force 
on the electron. This force is now along the direction of motion: 
the electron is accelerated. When the electron has crossed the 
interface into the Al side [d < 0), a negative charge is induced 
on the surface. Thus, a repulsive force Fjnd.s acts on the electron. 
This force is also along the direction of motion: the electron is 



accelerated. Therefore, the contribution of the induced surface 
charge to the force induced on an incoming projectile moving 
perpendicularly to the surface results in an acceleration at both 
sides of the interface. 

Regarding the effect of the induced surface charge on the 
probing projectile in the vacuum side, we can conclude the 
following: when a charged projectile approaches the inter- 
face to the medium perpendicularly, it is accelerated towards 
the surface. Conversely, when the projectile crosses the inter- 
face from the medium into the vacuum side, it is decelerated 
as it leaves the surface behind. This description is in accor- 
dance with Ref.'^^^ Regarding the medium side, it should be 
noted that there is an additional contribution to the 
induced force from the induced bulk charges. We anticipate 
that, as we will see below, the acceleration and deceleration 
of the projectile cannot be directly translated into energy 
gains or losses. 

In Fig. 5, we repeat the calculation for an oblique incoming 
trajectory with r/=120°. A general behavior is observed in 
psir, t) in the three studied cases {^ = 0, r/=120°, and r/=180°). 
As the projectile approaches the interface (regardless from 
which side), the magnitude of p^{r, t) increases monotonically. 
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Figure 5. Same as Fig. 3 for r;= 120° (oblique incoming trajectory). 



The projectile carries on with Its motion, crosses the 
surface and, as it leaves the Interface behind, the magnitude 
of Psix, t) does not decrease monotonically: It exhibits an 
oscillating behavior which is eventually damped out as the 
projectile moves further away from the interface. Conse- 
quently, the force induced on the projectile is also expected 
to oscillate as the projectile moves away from the surface. 
Such oscillatory behaviors will be seen repeatedly in the 
following sections. In the studied cases, the magnitude of 
Ps(r, f) Is substantially different from zero only within a very 
thin layer of a width In the order of 1 nm or 2 nm around 
the interface. Therefore, surface effects are to be expected 
only in a very thin layer around the interface of, at most, a 
few nm. This also Illustrates the surface sensitivity of con- 
temporary surface-analysis techniques: surface characteristics 
originating from a layer of a few Angstroems within the 
surface are resolvable. We point out how remarkable it Is 
that the method of image charges allows one to use the 
bulk dielectric functions of the media to describe processes 
that happen in a truly limited region of space. 



Stopping power 

In order to gain a better Insight Into the dynamics of a charged 
projectile in the vicinity of an interface. In this section, we examine 
in detail the variation of Its kinetic energy f per unit time: 



d£ d 1 , , dv ^ , . 

_ = __^(v.v) = mv.- = v.F,,(r,0 



Following Ref./^'' this equation can be trivially recast as 



d£ 

df ' 



V- Find(r, f) 



Zpedy^i^di 
dt 



Zoedy)i„d( 
dt 



(33) 



(34) 



which merely states that the variation In kinetic energy Is the 
variation In the total energy (kinetic plus potential, first two terms) 
minus the variation in potential energy (third term). On the other 
hand, the variation in kinetic energy of the projectile is given by 
the work W done on the projectile per unit time, which, in general, 
is the sum of a conservative and a dissipative contribution: 
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df 



df 



dt 



(35) 



The conservative contribution dW^ans^dt is given by the work done 
by the induced field, which is indeed conservative because it is the 
gradient of a potential: 



dWr, 



dt 



-Zne 



ind(r,f) 



df 



(36) 



This conservative work can be recovered in principle by returning 
the projectile to its original configuration. The remaining two terms 
in Eqn (34) constitute the dissipative part of the variation in kinetic 
energy per unit time: 



dl/Kdis 
dt 



v-Findfr,^ 



Zoed(/)i„d(r,f) 



df 



(37) 



This dissipated work is by definition not recoverable at later stages 
of the trajectory; it is invested in an irreversible way in energy losses 
(or gains in very exceptional circumstances, see below) of the 
projectile in the excitation of certain modes of the solid. Thus, it is 
precisely this part of the kinetic energy variation that is needed to 
sample genuine energy losses of the projectile in a Monte Carlo 
simulation. Consequently, in what follows, we will concentrate on 
the analysis of the dissipative part of the variation of the kinetic 
energy. Carrying out the full time derivative in Eqn (37), we obtain 



dWdi 



df 



^ v-Fi„d(r, f) 



-Zoe 



5Wnd("',f) 9<^ind('■,f)^>f/ 



af 



dx, dt 

where {xi, X2, X3} are {x, y, z}, respectively. Thus, we have 



J r=vt 

(38) 



dl^dis 
dt 



v-Fi„ri(r fl 



Zoe- 



dt 



V FinH(r,fl 



Zoe- 



,'9'/=ind(f,f) 



dt 



(39) 



The variation due to the acceleration or deceleration caused by 
Findlr f) (first term) is exactly cancelled out by the third term. 
Therefore, the conclusions drawn in the previous section regard- 
ing the acceleration or deceleration of the projectile should not 
be translated directly into an increase or a decrease in the kinetic 
energy, i.e. into non-recoverable energy losses or gains. It is the 
implicit time dependency of (p\ndif: f) which ultimately determines 
the dissipated kinetic energy. 

The approach we have followed (disregarding the conservative 
part of the variation in kinetic energy) is strictly valid for swift 
projectiles only. For these projectiles, the contribution of the 
conservative part along a reflected trajectory or a trajectory 
passing through a slab vanishes, provided that the velocity is 
not substantially modified during the interaction with the solid. 

Incidentally, we notice that in an infinite medium, there is 
space and time invariance and, thus, the induced potential 
evaluated at the position of the projectile cannot exhibit a 
dependency on time. Therefore, the dissipative part of the 
variation in the kinetic energy, Eqn (37), reduces to v- Find(i', f). 
Moreover, the conservative part of the variation in the kinetic 



energy, Eqn (36) vanishes. Under these special circumstances one 
can conclude whether the projectile undergoes an irreversible 
energy loss (or gain) alone from the direction of Find(r, f). In the 
case of a semi-infinite geometry, Eqn (39) must be used. 

The stopping power S is defined as the energy loss per unit 
path length. 



dWdis 
ds 



(40) 



where ds = vdf. This quantity has dimensions of force. Using 
Eqn (39), we have 



Zoe dipi 



ind 



r,f 



dt 



(41) 



In terms of the Fourier transform of the induced potential, ipi^d 
(q, (o)i we have 



.Zoe 1 f . r. 

S = I y / dq / d 



CO CO e''('' "-"Vind(q-«)- (42) 



Using Eqn (28), we have that the stopping power is given by 

[ 

dco CO / 

3 Jv? 



dq e''(''"-™)l^^(co-q.v) 



fa(q,co) f6(q,co) 



+ — d(co-qv) 
In 



pMd) ^ ps@{-d) 
£o(q,a)) Cb(q,(«), 



(43) 



Ps(q,«) 



0(d) @{-d) 



eo(q,w) £b(q,w) 



The first, the second, and the third summand are the stopping 
power of the projectile due to the interaction with the induced 
charges in the bulk, with the image charge (if selected), and 
with the induced surface charge, respectively. We will refer 
to these contributions using the symbols Sg, Sp, and Ss, respec- 
tively. In what follows, we shall concentrate on the surface 
model with (p,, p2, P3) = (0,1,0), which implies setting the sec- 
ond term of Eqn (43) to zero. This approach then encompasses 
the models of Refs.'^*'''^'^^'^^' The effect of the image charge, 
and to what extent it reproduces other models in the litera- 
ture, will be discussed in a future work. See also the note after 
Eqn (8). 

Let us first determine the stopping power due to the interaction 
with the induced bulk charges: 



2n^ V 

0(d) , 0{-d) 



dq / dco (u ^e'"''"-"''^(cu - q ■ v) 



(44) 



Assuming that £(q, a)) = £(q, co), which is typically the case in 
practice, the integral over the azimuthal angle of q can be trivially 
carried out to yield 271. The delta function can be used to carry 
out the integral over the cosine of the polar angle, c. A factor 
1 /qv appears from the change of variables in the argument of the 
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delta function and, since — 1 <(;<!, tine values of w range from 
— qv to qv. Finally, using the property £(q, a)) = £ (— q, — co), we 
obtain 



2{Zoef r.^r. 

Sg = — ^ — =^ / dq- dw oj Im 
7t Jo qJo 



(45) 



1 1 1 /■I 

Ss = — Uoef- / dq - dt) e"""!'^! / dos 
vJo ^ qJo 

X lm|e'*'^(ic?|||\/z| + q||-V||)f(q||,iq|||\/^| +q||-vj| 
0(rf) 0(-d) 



ea{q^q\\\vz\ + ^\\■y\\) £b(c?,iq||Kz| + q||-V||) 



(50) 



In an infinite medium, we would have only one of the dielectric- 
function terms, and we would retrieve the well-known expression 
for the stopping power in a bulk medium, see e.g. Eqn (1) of Ref '^^^ 
Similarly, for the stopping power due to the interaction with 
the Induced surface charge, we have 



S5 = ^(Zoe)^^ /■ dq \r do^co e^'C"-) 

271^ V Jul q^J-oo 



(46) 





r 0(d) 


0(-d)l 


c?2v|-h(cu-q||-V||)' 


.<:o(q.w) 


eb(q.«). 



Using the property e(q, tu) = f (— q, — a;), we can recast this 
expression as 



Ss = --{Zoer'-^ dq - 



dw 



X lm<^ e''*''"-™'^ q 



q|v2+ (^co- q||-V| 
0(d) 0(-d) 



(47) 



For times t < 0, we can carry out the contour Integration over co as 
described in the previous section. We define 53 = m — qn ■ V|| and 
recast Eqn (46) as 



I 1 \v.\ f ^ f°° CO + qii ■ Vii . , 

Ss = -^(Zoe)^^/ dq 4/ dS — ^e'^^V'-^ 



:f(q||,a)-Kq||-V||) 



0(d) 



0(-d) 



) -(< 



ta q, OJ + q,, • V|| fb q, « + qn ■ v 



We close the Integration contour through the upper half-plane In 
w and obtain 



Ss = f(Zoe)^l/dq l '^lll--l + ^ll-li e^.V^.M 



X f\^\\,'^q\\\vz\+%-^\ 

©(d) 



(49) 



0(-d) 



fia q,it?||Kz| + q||- V|| 1 eb q,iq|||\/z| + q||- V 



Using property 3, we obtain the final expression for the stopping 
power due to the Interaction of the charged projectile with the 
induced surface charge: 



for f < 0. For t > 0, we use Eqn (47). 

Figure 6 displays the sum of the stopping power due to the 
interaction with the induced bulk and surface charges, Eqns (45), 
(47), and (50), for an electron traveling near the vacuum-AI inter- 
face as a function of the distance to the interface. We consider a 
normal trajectory (left-hand panels) and an oblique trajectory 
(right-hand panels). In both cases, we consider the incoming and 
the outgoing trajectories. The stopping power for the incoming 
(outgoing) trajectory is displayed by the solid (dashed) lines. As 
we noted when discussing the induced surface charge, oscillatory 
behaviors appear when the electron has left the interface behind 
(both for Incoming and for outgoing trajectories). Thus, as Is well- 
known, there Is an asymmetry In the stopping power between 
the incoming and the outgoing trajectories: these differences 
become more important for trajectories crossing the surface with 
a direction close to the surface normal. 

The negative oscillation found in the stopping power for the 
outgoing trajectory in the vacuum side of the interface would 
strictly imply that the projectile gains energy in a certain range 
of distances at the vacuum side of the sample. However, the 
stopping power integrated over all positive distances* is positive 
and, thus, we can assume that the projectile always loses energy, 
thereby Introducing only minor errors. In the medium side of the 
interface, the stopping power Is always positive, and, therefore, 
the projectile is always slowed down In the medium. 

The magnitude of the stopping power at the vacuum side of 
the surface is generally smaller for the incoming trajectory than 
for the outgoing trajectory. A plausible explanation for this fact 
in connection with the conclusions of the previous section is 
the following. We recall that at the vacuum side of the incoming 
trajectory, the electron is accelerated towards the sample. Thus, 
the incoming electron spends less time In the vacuum side, 
and, therefore, less energy is lost In the incoming trajectory than 
in the outgoing trajectory. 



(48) Semiclassical approximation 



In the previous section, the stopping of the charged projectile 
has been attributed to the Interaction with the induced electric 
field, and therefore considered as a continuous energy-loss process. 
In reality, however, the stopping of the charged projectile occurs 
in the course of consecutive discrete Inelastic Interactions with 
the loosely bound electrons of the medium. The link between the 
classical and the quantized descriptions is made by means of the 
so-called semiclassical approximation, which consists In identifying 
the Fourier variables q and to with a discrete momentum transfer 
fi<{ and an energy loss fiw, respectively, in the course of an inelastic 
interaction. It Is now convenient to use the system of atomic units; 
h=e = m^=^, where m„ is the mass of the electron. Since In atomic 



*Note that this integration is only meaningful for very fast projectiles. 
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Figure 6. Sum of the classical stopping power due to the interaction with the induced bulk and surface charges (no image charges) for an electron near 
the vacuum-Al interface as a function of the distance to the interface, Eqns (45), (47), and (50). Here, d < 0 implies the Al side of the interface, whereas 
d > 0 implies the vacuum side. 



units, the energy loss Aoj reduces to oj, we use the same symbol co 
for the Fourier frequency and for the energy loss. 

The integration domains of q and co are now determined by 
the kinematics of inelastic collisions. Let the kinetic energy of 
the projectile be E before the collision and E—w after the 
collision. Let p and p , respectively, denote the wavevector of the 
charged projectile before and after the collision. The momentum 
transfer q is defined as 



P P 



Taking the square at each side of this expression, we obtain 



=p^ + p'^ - Ipp'cos e, 



(51) 



(52) 



where 0 is the polar scattering angle. In the non-relativistic 
case, energy and momentum are related by £ = p^/(2/D), namely 
p = \/2m£, where m is the mass of the projectile. Thus, we can re- 
cast Eqn (52) as 



= 2mE + 2m{E - W) - 4m ^E{E - oj)cose. 



If the projectile is an electron, we have (m= 1): 



q^ = 2E + 2{E -W)- 4v/£(E- W)cose. 



(53) 



(54) 



The maximum and minimum momentum transfer, q+ and q_, 
given by 



c?± = V2E±^2{E-W), 



(55) 



are achieved for 0 = 71 and 0 = 0, respectively. Thus, the integration 
domain of q will be restricted to [q^,q+]. This domain is valid for an 
electron both in an infinite and in a finite geometry: to derive it, we 
have merely used the definition of q and of the polar scattering 
angle. However, this does not imply that the azimuthal scattering 
angle is in general homogeneously distributed. Finally, the range of 
integration of the energy losses co must be restricted to [0, £]. 

Thus, in the semiclassical approximation, the stopping power due 
to the interaction with the charges induced in the bulk is given by 



2Z^ /'''+ 1 

Sg = J 6w CO I dq - Im 



0(d) , P2®{-d) 



£o(c?,C0) £b(a,w) 



(56) 



and the stopping power due to the interaction with the induced 
surface charge is given by 



Ss = --rZ^J-ii / 600 CO dq d&sM 

V Jo Jq_ Jo 



/•2n 

< / dip 
Jo 



qjv^+ I CO -q||-V| 



(57) 
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The Fourier phase might oscillate wildly and is therefore suscep- 
tible of posing numerical stability problems. In order to get rid of 
this problem, at least partially, we use lemma (Bl) to recast Ss as 



2n 



dip- 



E rq+ rn 

dw w dq d'd sini? 

0 Jq_ Jo 

1 



OJ — qn 



[0(-f)e"'"l''l + 0(f) (^2cos[t(^(u - q|| 
X lm<^ e'*°f q ' 



eo(q,a)) eb(q,a)) 



(58) 



DIIMFP 

In a quantum-mechanical calculation, one would derive a differ- 
ential interaction cross section for the inelastic scattering of the 
charged projectile with the electrons of the medium, doj/dcodq. 
The quantity (d(7/dojdq)dojdq is proportional to the probability 
that the charged projectile undergoes an inelastic interaction 
with an energy loss between to and co + dco and a momentum 
transfer within dq of q. The double differential inverse mean free 
path (DDIIMFP) is defined as 



d^i 
dojdq 



1 dui 
'yv'dwdq' 



(59) 



where Af is the average density of electrons in the medium. The 
quantity (d/(/dffiidq)dojdq is proportional to the probability that 
the charged projectile undergoes an inelastic interaction with 
an energy loss between co and w + dw and a momentum transfer 
within dq of q per unit path length. 

The electronic stopping power S is defined as the average 
energy loss per unit path length. To calculate it from the 
DDIIMFP, we would proceed as follows: 



dw CO dq 

0 Jq_ 



di?sini9 



dtp 



d/.i 
dcudq ' 



(60) 



where the angles j9 and tp are the polar and azimuthal angles of the 
momentum transfer in the Cartesian system of reference given in 
Fig. 1. Comparison of Eqn (60) with Eqns (56) and (58) allows us 
to identify 



d^i 



dcodq 



2n^ q 



Im 



0(d) , 0(-d) 



ea(q,co) £b(q,co) 



(61) 



and 



dcodq 



1 yl kzl 



1 



7l2 ° V 2 , 



X lm<^e'*''/^(q||,oj 



CO — q||-V|| 

0(d) 0(-d) 



Ca(q,co) C6(q,oj) 
[0(-f)e-'"l''l + ®(t){2cos\t[w - q||-V| 



(62) 



According to Eqn (53), the value of q is completely determined 
once a value for the energy loss co and the polar scattering angle 
0 is fixed. Hence, for a fixed primary energy, the DDIIMFP is a 
function of the energy loss and of the polar and the azimuthal 
scattering angles (the kinetic energy £, the depth d, and polar 
angle r] are parameters). The Cartesian components of the 
momentum transfer as a function of E, co, r;, and of the polar 
and azimuthal scattering angles 0 and 0 are given by:'^^' 



/ ±\/25 sin Tj ^ y/2{E — co)[cos ?/ sin 0 cosi^ + sin r/ cos 0] \ 



^yj2{E — co)sin 0 sinip 

\ COS f] + \/2(E - co)[sin r/ sin 0 cos i/j - cos cos 0 ] / 

(63) 



where the upper (lower) signs must be taken if r/ < 90° (77 > 90°). 

Figure 7 displays the contribution to the DDIIMFP due to the 
interaction with the induced surface charge [Eqns (62) using third 
component in Eqn (63)] as a function of the polar and the 
azimuthal scattering angles 0 and <^ for the following scenario: 
a 500 eV electron at a distance d= — 1 A at the Al side of the 
vacuum-AI interface moving with an off-normal angle 7/ = 60° 
(panel al) and 77=120° (panel bl). We consider a fixed energy 
loss CO =10.6 eV, corresponding to the excitation of a surface 
plasmon in Al. In this scenario, the DDIIMFP is significantly different 
from zero only for polar scattering angles which are less than ~ 2°. 
For higher energies, the DDIIMFP peaks at even smaller polar 
scattering angles. Hence, deflections in individual inelastic interac- 
tions can be neglected to a first approximation for £ = 500 eV. 
Furthermore, notice that the azimuthal scattering angle is not dis- 
tributed uniformly, as opposed to the case of an unlimited medium. 

In panels a2 and b2, we plot the absolute value of the 
momentum transfer in a direction perpendicular to the surface, 
\qX as a function of the polar and azimuthal scattering angles. 
There is a strong anticorrelation between the magnitude of jcjj 
and the intensity of the DDIIMFP: momentum transfers along 
the surface (scattering angles such that \q^ ^ 0) are strongly 
enhanced. Similar behaviors are found for other distances d at 
both sides of the interface, for other off-normal angles rj, and 
for other materials. 

Integration of the DDIIMFP over the kinematically allowed 
momentum transfers q gives a quantity that is proportional to 
the total distribution of energy losses in an inelastic interaction 
per unit path length. We will refer to this quantity as the DIIMFP 
and will use the symbol d/i/dco to refer to it. The area under the 
DIIMFP curve gives a quantity with units of inverse length, known 
as the inverse inelastic mean free path (IIMFP), Its inverse, the 
inelastic mean free path (IMFP), l-,, gives the average distance to 
the next inelastic interaction. 

From Eqn (61), we have that the contribution to the DIIMFP 
due to the interaction of the charged projectile with the induced 
bulk charges is given by 



d^iB(co) 
dco 



24 



dq -Im 



0(d) 



ea[q,o}) 



P2 



0(-d) 



fb(a,«) 



(64) 



These are the contributions to the DDIIMFP due to the interaction 
of the charged projectile with the induced charges in the bulk 
and at the surface, respectively. 



Similarly, from Eqn (62) we have that the contribution to the 
DIIMFP due to the interaction of the charged projectile with the 
induced surface charges is given by 
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Figure 7. (al and bl) Contribution to the DDIIMFP due to the Interaction with the Induced surface charge [Eqn (62) using third component in Eqn (63)] 
for a 500 eV electron at a distance d = — 1 A below the vacuum-AI interface for an energy loss oj= 10.6 eV as a function of the polar and azimuthal 
scattering angles, 9 and yn, for off-normal polar angles r; = 60° (al) and = 120° (bl). (a2 and b2) Absolute value of the momentum transfer perpendicular 
to the surface as a function of 9 and i/j for the scenarios corresponding to panels al and bl, respectively. 
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(65) 



Thus, we have expressed the DIIMFP of a charged projectile 
moving close to a planar interface in terms of the dielectric 
functions of the media at each side of the surface, of the velocity 
V of the projectile, of the distance to the interface, d, and of the 
polar angle of the surface-crossing direction, r;. 

All expressions derived In the present work require the bulk 
dielectric function of the medium as input information. The user 
can naturally select whichever dielectric-function model is more 
appropriate for his problem. In the results below, we have used 
a classical-oscillator-model dielectric function, using two different 
parametrizations. In the case of Al, we have used 



1 



f((7,w) 



qV4 



(66) 



where Zi is the atomic number of the target material, fj are the 
oscillator amplitudes (dimensionless), coy are the oscillator 
frequencies (in eV), and yj are the damping coefficients (in eV). 
The quantity is the squared plasma resonance frequency of 



the material. In the case of Al, we have used a single-oscillator 
model with a resonance at oji = 15.01 eV, amplitude fi=3, and 
damping coefficient =0.5 eV. For the 17 elemental metals listed 
in Ref.,'^^®' we have used the dielectric-function model described 
therein with the supplied parameter values. It should be noted that 
the fit parameter values given in this reference are those which give 
the best agreement between the computed DIIMFP and the DIIMFP 
retrieved from a deconvolution procedure from pairs of REELS spec- 
tra: a good agreement should be obtained between simulated and 
experimental REELS spectra, at least regarding the bulk losses. 
Finally, notice that other dielectric-function models can prove more 
realistic in practice, such as the Mermin dielectric function.'^'' The 
reader is referred to Ref.™ for a discussion of the effect of 
different dielectric-function models on the calculation of bulk 
energy-loss functions. 

Results 

Figure 8a displays the variation of the DIIMFP with d for a 500 eV 
electron moving perpendicularly from Al to vacuum (outgoing 
trajectory with 7/ = 0°), Eqns (64) and (65). When the electron is 
deep inside the Al side of the interface (solid black curve), the 
DIIMFP exhibits a peak at 15eV, the excitation energy of a bulk 
plasmon in Al. As the electron approaches the interface (solid 
red to solid blue curves), the bulk-plasmon peak is progressively 
suppressed while a peak is enhanced at 10.6eV, the excitation 
energy of a surface plasmon in Al. Thus, the effect of a planar 
interface is twofold: (i) it makes additional energy-loss modes 
available and (ii) it moves intensity from the bulk to the surface 
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Figure 8. DIIMFP for an electron near the Al-vacuum Interface, sum of Eqns (64) and (65). Examples for different kinetic energies E, depths d, and polar angles rj. 



energy-loss features. The latter effect is known in the literature as 
the Begrenzung effect (German for delimitation). In a future work, 
the separation of the surface DIIMFP Into a Begrenzung term and 
a clean surface term for the presented model will be studied in more 
detail. When the electron crosses the interface into the vacuum side 
(dashed curves in Fig. 8a), the bulk-plasmon peak is strongly sup- 
pressed, whereas the surface-loss peak prevails. This peak is in turn 
progressively suppressed as the electron moves further away from 
the interface (dashed blue through dashed black curves). The effect 
of the surface on the DIIMFP is noticeable in a region of about 10 A 
at either side of the interface. In Fig. 8b, the contents of Fig. 8a for 
the incoming trajectory (r/ = 180°) are reproduced. We see again 
that energy losses are possible at the vacuum side and that as the 
projectile moves deeper into the medium, surface-loss features 
are progressively suppressed to finally obtain bulk loss features. 

Let us consider the DIIMFP for the case d = 2 A (dashed purple 
curve), both in Fig. 8a (outgoing trajectory, interface already left 
behind) and 8b (Incoming trajectory, interface not crossed yet). 
The DIIMFP is larger in the case of the outgoing trajectory, for 



which the electron has already left the interface behind. This is 
also the case for other values of d at either side of the interface 
and conclude that the probability of undergoing a surface loss 
(which will be precisely defined below) is enhanced once the 
projectile has left the interface behind (at either side). 

Figure 8c and 8d illustrates the dependency on the polar angle 
of the trajectory, r/. We consider a lOOeV electron at cy= — 1 A 
(in Al) and d=^ A (in vacuum), respectively. The dependency 
on 7/ is not negligible: it is responsible for differences of up to 
100%, and It modifies the relative height of the bulk and surface 
peaks. The surface-loss peaks are Indeed enhanced for those 
values of r; corresponding to trajectories which have already 
crossed the interface: // > 90° for d = - 1 A and rj < 90° for c/= 1 A. 

Finally, Figs. 8e and 8f exhibits the dependency of the DIIMFP 
on £ for 7/ = 0° and d= — 1 A (Al side) and c/ = 1 A (vacuum side), 
respectively. In the Al side, the relative height of surface and bulk 
peaks varies significantly with E. In the vacuum side, we clearly 
see that the surface losses are enhanced for low kinetic energies. 
We can attribute this to the naive observation that for low 
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energies, the projectile spends more time crossing the surface- 
scattering region. Figure 9 displays a similar overview plot for the 
DIIMFP of an electron moving close to the vacuum-Au interface. 
The same conclusions apply. 

Comparison with literature models 

In Ref./^^' a similar calculation was carried out to describe energy 
losses of charged projectiles near a planar interface (we point out 
that the approach was rather different, in that a trajectory-dependent 
effective cross section was derived). The following simplification 
was made: 

T 0 obtain the last expression, it has been assumed that [. . .] 
£(q,co) = c(q^^,co) [. . .]. Although the validity of this approximation 
is not dear, it is needed to obtain an analytical expression.' 

This approximation (substitution of q bygy in the denominator of 
the expression for the induced surface charge density, Eqn (21)) 



will hereafter be called the Yubero et al. approximation. It was 
originally introduced for practical reasons in order to obtain tracta- 
ble expressions which could allow for an easy calculation of 
electron energy-loss properties near surfaces. To the best of our 
knowledge, no detailed discussion of the implications or the effect 
of this approximation exists in the literature. Moreover, several later 
models'^^'^^'"*^' use this approximation with the sole justification of 
a reference to the aforementioned work or to works which in turn 
make a reference to it. Another approximation that is made in some 
works (sometimes tacitly) consists in replacing q by q|| in the 
last bracket in Eqn (65). We shall refer to this approximation as 
the q qy approximation. 

Figure 10 shows the DIIMFP of a 500 eV electron moving 
perpendicularly from Al (d<0) to vacuum {d>0) for different 
distances d to the interface. The solid red curves are calculated 
using Eqn (65), whereas the dashed blue lines were calculated 
using the Yubero et al. approximation in the evaluation of the 
denominator of f(q||, w), and the thin solid black curves were 
calculated using both the Yubero et al. approximation and the 
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Figure 9. Same as Fig. 8 for Au and E=500 eV in the two upper rows. 
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Figure 10. DIIMFP of an outgoing 500 eV electron moving perpendicularly through the Al-vacuum interface evaluated at different depths d. The red 
curves have been calculated with the model described in this work, Eqns (64) and (65). The blue curves employ the Yubero ef al. approximation. The red 
curves, additionally, use the q qg approximation. The values of d are selected so that the thin black curve (Yubero ef al. and q qy approximations) 
can be compared with Fig. 2 of Reif.'^*' 



q qy approximation. The values of d were chosen so that the 
thin black curves can be compared with Fig. 2 of Ref.'^''' 

The effect of the Yubero ef a/, approximation is threefold: (1) it shifts 
the surface-loss peak to slightly smaller energy losses, ~ 0.6 eV in the 
considered case; (ii) it enhances the surface and the bulk losses by as 
much as 50% and 75%, respectively; (iii) it reduces the variance and 
the skewness of the surface-loss peak. The effect of the q^qy 
approximation is significant only for very close distances to the inter- 
face at the medium side. It enhances the surface losses (up to 
~20%) and suppresses the bulk losses (up to ~50%). Finally, in 
the vacuum side, the q^qy approximation has no effect, since 
e(q, oj) = 1 for any q. Notice that, in essence, both the Yubero ef al. 
approximation and the q^ q|| approximation imply replacing q by 
qy. We recall that, as we saw in the previous section, the momentum 
transfer in an inelastic collision takes place predominantly in 
directions along the surface. Thus, we conclude that both the Yubero 
ef al. and the q^qy approximations amount to an artificial 
enhancement of the surface-loss features. Finally, the difference 
between the DIIMFPs calculated with Eqn (65) and those calculated 
using the aforementioned approximations can be of the order of 
1 00% in the worst cases, e.g. the lower panels of Fig. 1 0. Notice that 
the differences in the position, shape, and width of the peaks in the 
DIIMFP will ultimately affect the quality of the agreement between 
simulated and experimental REELS spectra. The number of surface 
or bulk excitations, however, will depend only on the area under 
the peaks, which does vary (but not drastically) when using each of 
the examined approximations. 



In what follows, it will be shown how to derive the results of, e.g. 
Ref.'^^' from the general expressions given above. The contribution 
to the stopping power due to the interaction with the induced 
surface charge, Ss, will be compared with Eqns (10) and (11) in 
Ref.'^*"^ We consider Eqn (58) without restricting the integral to the 
kinematically allowed region: 
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Introducing the notation changes given in Appendix C, Eqn (67) 
becomes Eqn (C3), which can be recast more conveniently as 



cov. 



Ss = ^{Zoef- rdcof dQ „2 ^, , 



0(f)(2cos[5;f] -e-°W) 

£(Q,co) 1 -f(q,co) 
e(q, oj) 1 + £(Q, to) 



xlm 



-e-*l0(z) 



0(-z) 



e(z,Q,ffl) 



(68) 
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At tliis point, in order to reproduce Eqns (10) and (1 1) of Ref.'^**' 
we must assume that sign(t) = sign(z). This only holds for projectiles 
moving towards more positive values of z. Thus, we have 



1 1 f 

Ss = -^(Zoef- j dcoj dQ - 
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X Im 
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-Q|z|- 



^^ £(Q,co) 1 -£{q,co) 
£(q,aj) 1 +£(Q,co) 



(69) 



This expression matches Eqns (10) and (1 1) of Ref.'^'^^ 
Thus, the model in Ref.'^^*' has the following restrictions: 

• As pointed out by its authors, it is valid for outgoing projec- 
tiles, since we have had to impose sign(t) = sign(z). In recent 
works'^^'^^' a distinction is made between trajectories which 
move from the solid to the vacuum (labeled s^v) and from 
the vacuum to the solid (labeled v^s). 

• The integral over is carried out from — oo to +oo. This 
disregards the fact that there is a minimum and a maximum 
momentum transfer which are kinematically allowed in an 
inelastic interaction. As a result, arbitrarily large momentum 
transfers are included in the integral. Conservation of energy 
and momentum is therefore not completely satisfied. 

• In Fig. 10, in order to reproduce the figures in Ref.,'^''' both the 
Yubero et al. approximation and the q^qy approximation 
have had to be employed. 

Although this model is not the current state of the art, it is still 
occasionally used in recent works,'*^^ and it is simple enough so that 
it can be derived and examined within the present calculation in a 
straightfonward way. More detailed models exist, such as the model 
presented by Li et al. in Ref.,'^^^^ which makes use of both the Yubero 
ef al. and the q qy approximation and, moreover, it allows for a 
convenient separation of the surface contribution to the DIIMFP 
in a begrenzung term (suppression of bulk excitation modes) and 
a true surface-excitation term. A comparison between the model 
of Ref.'^^' and the model presented here is not straightforward 
and has been omitted for the sake of brevity: this comparison will 
be addressed in the forthcoming second part of this work. The 
performance of Li's model has been assessed by Novak, '"^^'"'''^ 
obtaining a remarkably good agreement between simulated and 
experimental REELS spectra. 

It should be emphasized again that the presented model 
comprises several models of the literature. By choosing appropriate 
values for (pi, p2, Ps), e.g. (0,1,0) for the Chen and Kwei model, and 
introducing further approximations as required, it is in principle 
possible to reproduce most of the subset of semiclassical models 
for surface excitations found in the literature.'^^'^^'"*^' As pointed 
above, the model of Yubero and Tougaard'^^' is not suited for 
comparison within the presented framework, for it takes a different 
approach. We do, however, point out that the Yubero and 
Tougaard model does account for the so-called interference effect 
between the outgoing electron and the field it set up upon 
entering the solid,'^^' whereas the model presented here does 
not. However, as shown by Vicanek,'^^' the interference effect 
might be significant for individual trajectories, but for a statistical 
ensemble of trajectories, for which the path-length distribution is 
usually much broader than the characteristic length for which 



interference effects play a role, the impact of these effects on the 
cumulative spectra is negligible. Last but by no means least, the 
agreement between the presented model and the Tung model'^^^ 
shall be investigated in the upcoming second part of this work. This 
model, albeit simple in that it is not position dependent, is neverthe- 
less of great use for practical applications in spectrum analysis. 

We point out that another class of models for the energy loss of 
charged projectiles near planar surfaces have been developed 
within a quantum many-body calculation such as those of 
Ding.'^^'^°' Their agreement with experimental data has been 
assessed,'^^' and their outcomes have been compared with 
those from models derived within the semiclassical dielectric 
formalism''**' such as the Li model'^^' exhibiting a very good 
agreement except at surface-crossing angles larger than 89°, implying 
that for most practical purposes, the semiclassical and the quantum 
computation are equivalent. These models have not been consid- 
ered in the present discussion, where we have concerned ourselves 
with derivations within the semiclassical dielectric formalism. 



Monte Carlo algorithm for the simulation of 
REELS spectra 

In this section, we describe a Monte Carlo algorithm for the 
simulation of the transport of charged projectiles in a semi-infinite 
geometry. See, e.g. Ref.'^^'"'''' for the sampling algorithm in the 
absence of surface excitations. In the present case, the procedure 
is very similar to the method described by Ding and Shimizu'"*^' 
and is based on the sampling algorithm of Coleman.'''^' It is 
assumed that a database of DIIMFPs for different materials, primary 
energies E, polar angles //, and depths d has been calculated on a 
grid of points of E, ;/, and d dense enough so that for a given mate- 
rial, a trilinear interpolation on E, q, and d is sufficient to interpolate 
values of the IMFP and the DIIMFP. Moreover, the differential cross 
section for elastic scattering (DCES), daJ6Q, is assumed to be 
calculated for a grid of energies £ which covers the range of interest 
for the simulation and which is dense enough so that linear (lin-lin 
or lin-log) interpolation suffices. The IMFP Ai{E,d,r]) and the elastic 
mean free paths l^iE) are tabulated beforehand. 
The core of the algorithm reads 

1 . Find the minimum mean free path among all processes (elastic 
scattering for all E, inelastic scattering for all £, d, rfi, a^^„. 

2. Initialize the trajectory (position and direction). In the case of 
(R)EELS, the primary electron is initialized at the highest posi- 
tion above the sample for which there are database entries. 
For XPS or AES, the electron is initialized at the corresponding 
emission depth. Notice that in the case of XPS, intrinsic losses 
due to the interaction of the photoelectron with the core- 
hole left behind are not accounted for by the present model. 

3. If the electron is in the surface-scattering region (between the 
minimum and the maximum d for which there are database 
entries): 

(a) Take a steps ~ (1 /Amin)e"^^''"™ along the current direction. 

(b) Interpolate and at the current kinetic energy, position, 
and direction of flight (/le = 0 at the vacuum side). 

(c) Sample an event of the type: 
• null with probability A^Jn/(^.min + ^= ^ + ^-^^ V 



• elastic interaction with probability a^ /{a 



\. Let £,<£<£,+i. Sample a deflection from d(Te(£/)/dQ 
or df7e(£;+i )/dO with probabilities 
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and 



ln(£,-+i/£) 
ln(f,+i/£/) 



ln(£/£/) 
ln(£/+i/E,) 



(70) 



(71) 



respectively.'"*^' 
ii. Update the direction of flight of the 



electron. 



inelastic interaction with probability /IrVl'^min + + 
1^^) ■ We now need to sample an energy loss from 
the database for the current E, ij, and d. Let Ei<E<Ej+t, 
Tjj^riKrij+i and dj^<d <dj^+^. LetXdenote £, r/, or d and 
let the subscript a denote the corresponding subscript 
/, j, or k, as needed. We take X=X^ or X=X„ + i with 
probabilities 



and 



Xa+1 — X 



X — Xa 
Xs!+1 — Xc 



(72) 



(73) 



respectively. Let a,-, 09, and denote the chosen value for 
the subscripts and k, respectively. Once the active sub- 
scripts a,-, and are determined, we sample an energy 
%.,d^i^ ) and update the energy of the 



loss from n\^co;Ert- 
projectile. 

Else, if the projectile is outside of the surface-scattering 

region: 

(a) Interpolate A; and 1^ at the current energy. Calculate the 
total mean free path a^'^ = a^^ + 

(b) Take a step s ~ (1 along the current direction. 

(c) Sample an event of the type: 

• elastic interaction with probability /i,gV(''-e^ +''T^)- 

i. Let E;<£< £,+1. Sample a deflection from d(Te(£/)/dn or 
d(7e(£/+i)/dn with probabilities k, and 71,+ , defined 
above, respectively. 

ii. Update the direction of flight of the electron. 

• inelastic interaction with probability X^^ /[a^^ + /'-i"^). 

i. Let Ei<E< £;+i. Sample an energy loss from n{E^ or A((£;+i) 
with probabilities 71, and 71,+ 1 defined above, respectively. 

ii. Update the energy of the electron. 

4. If the electron leaves the sample (surface-scattering region) 
without entering the analyzer or if the electron leaves the 
energy window of interest, disregard the trajectory and 
sample the next one. Otherwise, sample the current trajec- 
tory further by going back to step 3. If the electron leaves 
the sample in a direction within the solid angle of detection 
around the detector direction, update the corresponding 
histograms. 

Comparison with experimental REELS spectra 

The algorithm outlined in the previous section has been used 
to simulate a database of experimental REELS spectra'^''^"^^' for 
24 elemental solids measured on a VG Microlab UHV system 
equipped with a hemispherical analyzer, used with a constant 
resolution mode with a pass energy of 20 eV. The beam of 



incident electrons impinged normally on the sample and the 
analyzer formed a polar angle of 60° with the surface normal. 
The pressure during the measurements was in the 10^® mbar 
range. See Refs'"*^"^*' for more details. The experimental spectra 
were converted to absolute units of 1/eV through a division by 
the area under the elastic peak and by the bin width. 

The database of DIIMFPs and DCESs required for the simulation 
has been generated as follows. The optical data provided in 
Ref.'^^' has been used to generate the required database of 
£-, r]-, and d-dependent DIIMFPs, using Eqns (64) and (65). We 
have used the following grid: £={100,500,1000,5000} eV, 
r/ = {0°,60°,120°,180°}, d = {±1 5,± 5, ± 2,± 1,± 0.5,± 0.1} A. Addi- 
tionally, a database of 50 bulk DIIMFPs has been generated for 
energies equally spaced on a logarithmic grid from £=100 eV to 
£=10 keV. DCESs for 150 energies equally spaced on a logarithmic 
scale from £ = 1 00 eV to £ = 1 0 keV have been generated with the 
ELSEPA code,'^^^' using a muffin-tin model potential obtained 
from the Dirac-Fock electron density, with exchange effects 
accounted for by means of the Furness-McCarthy exchange 
potential. Correlation-polarization and inelastic-absorption cor- 
rections have not been considered. 

The simulation of the aforementioned database of experimental 
REELS spectra (as well as the calculation of the DIIMFP database) 
was carried out at the Vienna Scientific Cluster, using Intel Xeon 
X5550 (Nehalem) processors with a 2.97 GHz clockspeed. The 
sampling algorithm attained an average speed of 5 x lO"* up to 
10^ simulated trajectories per second; 10^ trajectories were 
sampled for each REELS spectrum with a resolution of 0.1 eV. In 
Fig. 11, we compare six of the simulated REELS spectra with their 
experimental counterparts, with primary energies from 700 eV to 
1 200 eV and metals ranging from Ti (Z= 22) to Bi (Z= 83). The solid 
curves are the simulated spectra, whereas the dashed curves are 
the measured spectra. Notice that the comparison is in absolute 
units of reciprocal eV: the experimental spectra have been divided 
by the area of the elastic peak (which was subsequently subtracted) 
and by the bin width; the simulated spectra have been divided by 
the number of elastically detected electrons and by the bin width. 
The simulated spectra are in good agreement with the experimen- 
tal spectra, the discrepancies being between 5% and 10% in the 
considered energy-loss range. The agreement is strongly depen- 
dent on the quality of the used optical data, to a much higher 
degree than in the simulation of bulk losses only. Finally, we point 
out that the intensity scale of the normalized experimental spectra 
might be off by ~ 5% due to the uncertainty in the determination 
of the area under the elastic peak. We point out that, although 
the primary energies in Fig. 11 vary between 700 eV and 1200eV, 
a similarly good agreement was found for other energies. 

As outlined above, the algorithm requires considerable 
computational effort. If only bulk losses are considered, the 
trajectory reversal algorithm'^^' can be employed to drastically 
decrease the simulation time. However, this algorithm relies on 
the symmetry of the loss characteristics between the incoming 
and the outgoing trajectory. This symmetry is lost in the 
surface-excitation model described here. The trajectory-reversal 
algorithm can therefore not be employed to accelerate the 
simulation. In a future work, it will be discussed to what extent 
the present model can be used to derive a cumulative surface- 
loss distribution (the so-called differential surface excitation 
probability (DSEP)), an approach similar to that of Ref.'^^' DSEP-type 
loss distributions can be easily embedded in the partial-intensity 
approach'^"" for a simplified and yet accurate description of surface 
excitations in electron spectroscopy. 
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Finally, in Fig. 12, the effect of the Yubero et al. and the 
q-^ qil approximation on simulated REELS spectra is exposed. We 
consider the REELS spectrum of 500 eV electrons impinging 
perpendicularly on an Al surface with a ring detector accepting 
outgoing trajectories with off-normal polar angles between 
58° and 62°. The solid curve was calculated using the model 
described in this work, whereas the dashed curve was calculated 
using the model presented in Ref.'^*"' (the model presented 
here was used with the approximations addressed above to match 
the model in Ref.'^*"' The database of calculated DIIMFPs comprised 
the following parameter values: £ = {450,500}, r/ = {0,60, 120, 180}, 
c/={±15,±5,±2,± 1,±0.5,±0.1}. The differences we observed 
in the DIIMFPs between the two models are naturally translated 
into the REELS spectrum: surface features are accentuated and 
shifted towards smaller energy losses, and the relative intensity of 
bulk and surface losses is modified. An experimental REELS spec- 
trum of 500 eV impinging normally on Al and being detected along 
a polar angle of 60 ± 1 2° is shown. The relative bulk/surface height 
and the asymmetry and width of the surface-loss peaks are well 
described by the present model. Regarding the intensity scale 
and the overall sharpness of the simulation results, it must be 



pointed out that the simulation of REELS for Al is typically problem- 
atic with the employed dielectric-function model: the DIIMFPs 
derived from optical data usually give too narrow loss features in 
comparison with the experimental loss features. Similar difficulties 
are found for other nearly free-electron metals such as Si. Neverthe- 
less, they are ideal candidates for a benchmark of model calcula- 
tions due to their sharp loss features, both in the bulk and in the 
surface. A convolution with a normalized Gaussian peak of the 
estimated experimental energy resolution might lead to a change 
in the relative intensity of surface and bulk peaks. We note that 
the employed dielectric-function model does not properly include 
electron-hole pair generation. The use of other dielectric-function 
models which account for electron-hole pair generation such as 
the Mermin dielectric function'^^'^°' should lead to a better agree- 
ment between simulated and experimental REELS spectra for 
nearly free-electron materials such as Al and Si. 

Conclusions 

A comprehensive and self-contained description of the dynamics 
of non-relativistic charged projectiles in the vicinity of a planar 
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AL E = 500 eV, 9^^ = 0°, 6lout = G0° ± 2°, ring detector 

0.9 I , , , , , , 1 0.15 

j : ; : Present model — 




Energy loss (eV) 



Figure 12. Simulated REELS spectrum of 500 eV electrons impinging perpendicularly on a polycrystalline Al sample detected by a ring detector at 
60° ± 2° using the model described in this work (solid curve) and the model in Ref." (dashed curve, see text). The dot-dashed curve corresponds to 
an experimental REELS spectrum measured for the same geometry (with a HMA at 60° off-normal angle with a half-opening angle of ~12°) and 
primary energy. 



interface is provided within the dielectric formalism, considering 
only the longitudinal part of the force induced on the projectile. 
The presented framework allows one to derive a differential IIMFP 
(DIIMFP) of a charged projectile as a function of its speed, of its 
direction of motion, and of its distance to the interface. The deriva- 
tion has been carried out in such a way that (i) it encompasses a 
number of models in the literature and (ii) unnecessary simplifica- 
tions are avoided. The effect of the assumptions and simplifications 
made in literature models has been investigated, revealing differ- 
ences of up to 100% (in the worst cases) in the intensity of the 
DIIMFP, as well as appreciable differences in its shape. DIIMFPs 
calculated with the model described here have been used as input 
for a Monte Carlo algorithm for the simulation of REELS spectra, 
obtaining good agreement between simulated and experimental 
spectra in absolute units. We therefore conclude that the presented 
model can be applied in a detailed simulation of energy-loss 
processes of charged projectiles near planar interfaces. 
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Lemma 1 

If f(co) is a function of tu with no poles in the upper complex half- 
plane, then 



doj e 

, f{co) I J-^ °'- + '"^ 



for t<0, 



dw [2cos(cof) - e-l'l°l for t > 0. 

(B1) 

For f<0, the lemma follows from direct integration by residues of the 
integrands at both sides of the equality along a semi-circular path 
on the upper half-plane, picking up the pole at oj = '\a. For f>0, 
we can use the trivial identity 



2cos(ta;) - e'^ 



(62) 



and close the integration contour again along a semicircle on the 
upper half-plane. 



Appendix A: Fourier transform 

We shall use the following definition of the direct and inverse 
Fourier transform and its inverse, using symmetrical factors of 
1 /V2n for each of the four integration variables: 

f{q,aj)=-\f dr Tdt e-'C'-'^'VlrO (A1) 

(27r) Jk' J-oo 



1. Properties of the dielectric function for complex arguments 

When carrying out contour integrals over oj, it becomes necessary 
to evaluate c{q,oj) and fiq^.w) for complex values of oj. The proper- 
ties of these functions for complex oj will now be examined. 
Since in all practical cases, we consider c{q,oj) = e{q,oj), we shall 
restrict ourselves to a classical-oscillator dielectric-function model 
of the form of Eqn (66). 



f{r,t)=-^ f dqT dw e'C-"') f (q, w) . (A2) 

2?! J«> J-oo 



As an example, consider the Fourier transform of a charge 
distribution p{r,t) = Zoe8{r — vf): 



p(q,co)=-^/ dr/ df e-'(''--"<5(r-vt) 

{2k) Jw' J-oo 



7t) J-oo 



(2 
Zoe 

" 271 



(A3) 



where we have used a familiar representation of the delta 
distribution. 



Appendix B: Mathematical tools 

To evaluate integrals which contain a phase susceptible of 
oscillating wildly, we propose the following pre-conditioning 
process. 



Property 1 

For a complex frequency of the form co = i(7|||i/J -Fq||V||, we have 



1 _ 1 

q||V|| + i£?||kz|) e*(^c?, -q||V|| + jc?|||i/z| 

Proof: For (M = q||V|| + ic?|||vj, we have 
1 

1 

e(^q,q||V|| + iqiil^z'' 



(B3) 



^«/ + '7V4-((q||V||)'-q|vl) 



i)Vq||V|| +)';q|||Vz| 



(B4) 



If we change qy by — qy, we have 
1 

^= 1 



f q,-q||V|| + iq|||i'z| 



7" tu/ + <tlA - (^(q||V||)^ - q\vlj + i-/yq||V|| + yfl^\\v,\ 



1 

f*('7>qi|V|| + \q\\\vz 
Property 2 



(B5) 



□ 



nqil.iQllkzl + q||V|| 





P2+P3 


1+Pl 


fb(o,i£?|||i'z|+q 


V||) £„(0J()|||\/,|+q||V||) 


roc 

1 dKz -i:!-! 


1 1 1 


Co(<J||.Kz. 


flllkzl+qilVll)) ' Ci,(lJ||,Kz,ilJ|||l/z|+q||V||)) 



(B6) 
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Proof: From the definition of f{q\\,oj), we have that 

o) — q||V|| 



(B7) 
□ 



For CO = iq|||vJ + q||V||, we have that /Cz= i<7||Sign(vJ. Therefore, 
c/2 = q| + /c| = 0. 

Property 3 

For a complex frequency of the form tt) = i(7|||i/J + q||V||, we have 



{B8) 



M 111 ' lll^ll + '^lll^'^l )=''*( "111 ' "lll^ll + '"^ll l*'^ 



□ 



provided that e{q,co) = c{q,oj). 

Proof: It follows directly from Eqn (B3) and the definition of f(q\\,w). 

Appendix C: Notation changes in the literature 

We introduce the following notation changes in order to comply 
with Ref.: "'^i 



P2 = T Pi = P3 = 0, eo(q, oj) = T <:b(q, 03) = £{q, a), 
Q = q|h ffl = a)-q||V||, d = z, 



(u-q||V|| 



1 



1 1 



(CI) 



e(Q, w) K ^ e(q, oj) 



e(z, 0,0)) 71 



1 1 

q2 w) 



With this new notation, we have 



^(q||,«)=- 



1 



Q e{Q,oj) 1 -£(q,m) 
71 e(q,tu) 1 +e(Q, co) 



{C2) 



Thus, 



^(Zoe)'l / do) / dQ ^ 



Co|Vz| 



0(-t)e-°l^l + 0(f)(2cos[ajr] - e 



x^lmHS 



Q,tQ) 1 -£(q,m) Z'" 
q,oj) 1 +f(Q,oj) ■_/__, 



d<7z 



0(z) 



e(q.™) 



{C3) 



00 
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